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Abstract 

We have employed the steepest descent method to optimise the variational ground state quan- 
tum Monte Carlo wave function for He, Li, Be, B and C atoms. We have used both the direct 
energy minimisation and the variance minimisation approaches. Our calculations show that in 
spite of receiving insufficient attention, the steepest descent method can successfully minimise the 
wave function. All the derivatives of the trial wave function respect to spatial coordinates and 
variational parameters have been computed analytically. Our ground state energies are in a very 
good agreement with those obtained with diffusion quantum Monte Carlo method (DMC) and the 
exact results. 



1 



I. INTRODUCTION 



Quantum Monte Carlo (QMC) method has constituted an efficient and powerful nu- 
merical method for solving time-independent many-body Schrodinger equation mainly in 
chemistry and solid state physicai'^'^'^i^'^'^i^i^i'i^. Among various approaches to QMC namely, 
random walk, diffusion. Green-function etc, variational quantum Monte Carlo (VMC) has 
been extensively studied in recent years-. In VMC method, a parameterized many-body 
trial wave function is optimised according to Raleigh-Ritz variation principle. In practice, 
this task is done utilizing a numerical algorithm for optimisation the parameters. Vari- 
ous algorithms have been proposed and implemented in the framework of QMC such as 
Newtoniiii^ii^, steepest descent (SD)^ii^ii^, perturbative optimisationi^ii^ and linear opti- 
mization method^ii^'22. The wave-function optimisation is implemented via two schemes 
namely energy minimisation and variance minimisation. These methods have their own 
merits and disadvantages. A basic task in VMC is the evaluation of first and second deriva- 
tives of the local energy = ^ respect to variational parameters and spatial coordinates 
or a combination of them is the trial wave function). Despite normally the first deriva- 
tive is analytically evaluated and the second derivatives are calculated numerically^ there 
are papers in which second derivatives are also calculated analytically^. Numerical evalu- 
ation of second derivatives causes a systematic error into the problem. To the best of our 
knowledge, the SD method has only been utilized in the variance minimisation approach^^. 
Our objective in this paper is to show that implementation of the SD method in the direct 
approach of energy minmisation yields reasonable results. We report our results for the 
ground state energies of He, Li, Be, B and C atoms and compare them to the results in the 
literature obtained by other methods. 

II. VARIATIONAL WAVE FUNCTION AND STEEPEST DESCENT OPTIMI- 
SATION METHOD 

A. Theoretical background 

Let us briefly explain the basic ingredients of the VMC method. In the VMC method 
a trial many-body wave function \l/(-R, {cm}) containing a set of M variational parameters 
Ci, C2, ■ ■ ■ , Cm is considered. R denotes the position set of electrons. We confine ourselves 
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to Born-Openheimer approximation in which the nuclei are assumed static and only the 
electronic degrees of freedom are taken into account. The parameters Cm are varied according 
to the Raleigh-Ritz variation procedure so as to minimize the variational energy E{{cm}) 
defined as follows: 

j^.r ^ J^*{R,{c^})H^{R,{cm})dR 
J^*iR,{cm})^iR,{c^})dR 
Where H is the many body system Hamiltonian. We ignore relativistic correction and 

take the Hamiltonian as follows (in Hartree atomic units): 

^ = -^Ev.^-E^ + E;^ (2) 

i i,I *^ i<j *i 

Small letters refer to electrons and capital ones to nuclei. Zj is the electric charge of 
the I-th nucleus and Vij denotes the distance between electron i and electron j whereas 
Tj/ denotes the distance between electron i and nucleus /. Moreover, we restrict ourselves 
to real-valued wave function \1/ and omit the complex conjugate symbol afterwards. By 
introducing a local energy = ^ and a normalized probability distribution function 



p{R^ {cm}) = (Z2^^'\'^"'\\j^ We recast equation (1) in the following form: 



E{{c^]) = j p{R,{cm})ELdR (3) 
It is now possible to approximate the integral by the standard Monte Carlo procedure: 

Ni^iC 

E{{cm}) = lirriN^.c^oo^ — ^ El^^ (4) 

iVAfC s=l 

El,s denotes the local energy of the s-th sample of the configuration-space and Nmc is 
the number of Monte Carlo sampling for evaluation of the integral. 



B. Optimisation of the wave function: steepest descent method 

The next step is to find the optimal values of the parameters which minimise the objective 
function i.e.; the variational energy EM-'^'^. There exists numerous optimisation methods 
such as Newton^i'^'^, steepest descenlr^'^'^, conjugate gradient etc in the literature. Here 
we focus on the simplest of them and show that despite simplicity algorithm is capable of 
exhibiting a satisfactory performance despite not receiving much attention. We briefly recall 
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the main ingredient of this method. Having numerically computed the energy E for a given 
set of parameters in (4), we iteratively update the values of the parameters according to the 
following procedure: 



The vector c'l" = (ci, C2, ■ ■ ■ , cm) denotes the parameters, k is the iteration step and a 
denotes the constant of the SD method. The vector g"!" = (|^) f^j ■ ' ' i ^fp) gradient 
vector of energy respect to the parameters. We note that in some cases we should vary the 
SD constant in each iteration step to get the desired optimum value. In order to utilize SD 
method, we should evaluate the energy gradient vector. This has been done in details 
We only quote the result: 

hrriNMc^oo^T — 2^ {[■^r—)s[EL,s - E)) (6) 



dc 



In eq. (6) (^^)s denotes the logarithmic derivative of wave function evaluated in the 



s-th MC configuration. If the constant a is appropriately chosen the sequence converges 
to c* after some iterations. 



C. Trial wave function and its parameters 

We wish now to introduce the structure of the trial ground state wave function we have 
implemented in our calculations for simple atoms. We have taken the following well-known 
form for \EtS^i^: 



^ = D^D^e-^ (7) 

In which D"^ and are up-spin and down-spin Slater determinants and J is the Jastrow 
factor. The number of spatial orbitals Nup{Ndown) in the construction of Slater determinant 
D^{D^) equals the number of spin up (down) electrons and depends on the atom we consider. 
Note that Nup + Ndown = N where is the number of electrons in the atom. For the basis set 
in the construction of up and down Slater determinants we have used a variant of Slater-type 
s and p orbital as follows^^iS^: 



4 



(r) = E Cke-<^'' + E Cire-'^'^' (8) 



k=l k=l 



,r 



E xDke'^^'^ (9) 



k=l 



Analogous definitions goes for py and orbitals. We liave set P = 3 in all our calcu- 
lations. Henceforth, the parameters are Ck,C'i^,Dk,Ck,Ck where k = 1,2,3. For the 
Jastrow factor we have taken the following form: 

J = E (10) 

i<j 

The sum goes over all the particles (electrons) and Uij has the following dependence on 
distances: 

Tj is the distance between electron i and the nucleus, r^j is the distance between electrons 
i and j. Exponents m, n, o are positive integers and the sum over mno denotes the sum over 
given values of these integers. We adopt the following choice of integers^^ m, n, o: 

{(0, 0, 1), (0, 0, 2), (0, 0, 3), (0, 0, 4), (2, 0, 0), (3, 0, 0), 

(4,0,0), (2,2,2), (2,0,2)} (12) 

Equation (11) includes electron-electron correlations (terms with m = n = 0), electron- 
nucleus correlations (o = as well as one of m or n zero) and also electron-electron-nucleus 
correlations ((2,0,2) and (2,2,2) terms). Here we have considered the simplest choices com- 
patible with electron-electron and electron-nucleus cusp conditions. The origin of three body 
correlation terms in (11) stems in the back flow correlation firstly suggested by Feynman 
and Cohen^. We refer the readers for more details to reference [25]. The Jastrow function 
has nine independent parameters Cmno- Each s type orbital contains twelve parameters and 
in a p orbital we have six parameters. We note that after imposing electron-nucleus cusp 
conditions, one parameter from each s orbital will be fixed. 
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D. Variance minimisation method 



In the preceding sections, we outlined the basics of energy minimisation method. In this 
method, we minimise the variational energy. In recent years, an alternative scheme the 
so-called variance minimisation has been introducecy^i^^ and has become one of the most 
frequently used method in the literature. This method has shown to provide some advantages 
over the straightforward energy minimisation. We now briefly review this method. Instead 
of energy, we minimise the variance of the local energy El—^'- 

/*2(B,{c„})<ifi 

All the other steps are analogous to those in the energy minimisation method. To imple- 
ment the SD procedure we only should replace the energy gradient vector with the variance 
gradient vector. Derivatives of o"^ respect to parameters have been evaluated Here for 
simplicity we use the following expression which ignores the change of the wave function^: 

(14) 

More concisely the above approximation corresponds to underweighted variance minimi- 
sation method. The average is taken with the normalised probability function p{R, {cm}) = 
J^^{^'{c"^})dR ' approximate the average in (14) by a sum in MC approach. Note 

that when implementing this method, we have to replace E with in the gradient vector 
g"!" in equation (5). In the next section our results will be reported. All the computational 
details of the calculations are explained in the appendix. 



III. APPLICATION TO ATOMS AND DISCUSSION 

We have implemented the steepest descent optimising method to find the ground state 
energy and wave function of atoms He, Li, Be, B and C by two approaches of energy and 
variance minimisation. Let us now explain our procedure of energy minimisation. It consists 
of three steps: anticipating the variational parameters, finding the optimised value of steepest 
descent parameter a and eventually the fine tuning of variational energy. Step one begins 
with random initialization of the parameters values. Initial values of Jastrow parameters 
are randomly chosen in vicinity of zero. We then set the SD constant a to a rather high 
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value say a — 0.1. Next we proceed with some iterations of (5) until the variational energy 
reaches approximately to the exact ground state energy. This comprises step one. During 
this step integrals in (1) are evaluated by the standard MC Metropohs method. Each MC 
move consists of a random selection of an electron and displace it from its position by the 
vector d. The move size which is the length of 6 is randomly chosen (uniformly) from the 
interval [0,^r„aa;]. The direction of 5 is uniformly chosen between zero and 2tt. We took the 
number of Monte Carlo steps Nmc equal to 3 x 10^. We discard the first 10000 steps to 
ensure reaching equihbrium. Averages are separated by 20 MC steps to suppress the effects 
of correlations among generated MC configurations. The MC maximum move size 5max has 
been typically 0.3 (Hartree atomic units) with the acceptance ratio around 70 percent. At 
the end of step one, which normally takes 5 — 6 iterations, variational parameters should 
have reached to the vicinity of their ultimate values. We put their latest values in the code 
and re run it. This is the beginning of step two. We then proceed with some iterations 
until the iteration series of the variational energy begins to diverge. This shows that by the 
current value of the SD parameter we can no more reach the true energy. Here we reduce 
a to a smaller value say one order of magnitude smaller and repeat the procedure until the 
iteration series of energy begins to diverge or strongly oscillates. Wc repeat this a reduction 
procedure until further reduction of the SD constant a does not lead to divergence of energy 
iteration series. Normally after 4 — 5 repetitions we achieve our aim and a reaches to a value 
of the order 10~^. This marks the end of step two and by now we have an iteration energy 
series. In figure (1) we have depicted such series of Be ground state energy obtained in the 
method explained above. Corresponding series for other atoms are similar in nature. 

It is seen that after roughly 20 iterations we reach a steady state regime. Next in figure 
(2) we exhibit the dependence of absolute value of VcE on the iteration number. 

We see that iVc-E] tends to a small value for a sufficient number of iterations. Theoreti- 
cally it should goes to zero but due to numeric computations it does not approach to zero. 
Finally in third step, in which the fine tuning of energy in performed, we choose the varia- 
tional parameters obtained from that iteration in the series which has the smallest energy 
and re run the code for these values of parameters for a longer MC run of Nmc = 2 x 10^ 
steps to find the energy in a fine tune manner. All the reported vahies in table I have been 
obtained in this way. In figure (3) the energy error is shown for both methods of energy 
and variance optimisation. We recall that the error has been obtained from the following 
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FIG. 1: Ground state energy of Be in two methods of energy and variance minimisation. 




FIG. 2: |Vc-E| vs iteration number. 

formula: 



error 



\ Nmc - 1 



(15) 



We see that the error is lower in variance minimisation method. In both methods, the 
error decreases as we increase the iteration number. In table I, we report the ground state 
energies we have obtained and compare them both to the existing computational results 
in the literature obtained by other methods^'^^'^ and the exact ones'^"^. The comparison 
shows that the steepest descent method is capable of minimising the energy to a very good 
precision. In fact our results by energy minimisation method is in most of the cases even 
better than those reported Besides Be, in all the atoms, our energy is comparable 
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FIG. 3: QMC error vs iteration number. 

to the DMC energy. This marks the efficiency of the steepest descent method at least for 
hght atoms. We should like to emphasize that simplicity is the main merit of our approach 
which can turn it into an efficient method al least for simple atoms. Our results obtained by 
variance minimisation method is less favourable in comparison to our energy minimisation 
ones, however, the error in the energy minimisation method is larger. Comparison of our 
variance minimisation results to those in-- (which have been obtained by VMC in Newton 
optimisation method ) shows that the results of^-*^ is slightly better than ours. We have also 
compared our results to the recent paper of Brown et aF^ in which besides single determinant, 
multideterminant trial wave function have been employed. It is seen that the accuracy of 
our results is comparable to those exhibited in^^. 

IV. SUMMARY AND CONCLUDING REMARKS 

In summary we have applied the steepest descent optimisation method to optimise the 
parameters of the QMC many-body wave-function in some light atoms. Two schemes of 
energy and variance minimisation have been implemented. The key features to achieve the 
correct minimum is to vary the SD constant a appropriately. Our results are in a well 
agreement with exact results and those obtained by DMC. We note that all the derivatives 
of the trial wave function respect to spatial coordinates and variational parameters have 
been analytically calculated. 
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TABLE I: Variational energy with error bar for atoms He to C (All the energies are in Hartree). 





He 


Li 


Be 


B 


C 


i?(energy minimisation) 


-2.9037(4) 


-7.4780(4) 


-14.648(1) 


-24.640(9) 


-37.831(8) 


-^(variance minimisation) 


-2.9031(2) 


-7.4757(3) 


-14.6443(9) 


-24.6244(3) 


-37.807(5) 


EoiRem) 


-2.903719 


-7.47806 


-14.66736 


-24.65391 


-37.8450 




-2.903717(8) 


-7.47722(4) 


-14.6475(1) 


-24.6257(1) 


-37.8116(2) 




-2.9029(1) 


-7.4731(6) 


-14.6332(8) 


-24.6113(8) 


-37.7956(7) 


EDMciRem.) 


-2.903719 


-7.4780(2) 


-14.6565(4) 


-24.63855(5) 


-37.8296(8) 


EvMciRem) 


no report 


-7.47683(3) 


-14.6311(1) 


-24.6056(2) 


-37.8147(1) 
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VI. APPENDIX 

In this section we give some details of the manipulations for the evaluation of the integrals 
(3) and (14). In evaluation of these integrals, we have analytically calculated two basic 
quantities El, and Let us first consider El. According to its definition we have: 

The first term is kinetic energy and V represents the potential energy. Calculating V is 
straightforward. To calculate the kinetic energy KE we rewrite it in the following form^: 

^^ = E-^[v?in^ + (^^i^^)'] (17) 

Concerning the form of the trail wave function, \1/ = D'^D^e'^ we have: 

In v]/ = -^^^D'^ + -^^rD^ + J (18) 



and 
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l^VlD^ + Vp (19) 

We note that the entry Dij of any determinant D equals (pi{rj) in which the 
orbital (pi is the 2th orbital in the construction of D. To evaluate the matrix 
"^.[D wc only have to replace the ith column of matrix D Cj by a new column 

ND{^Nn)V- The matrix Vf-D is analogously constructed 
but with Ci = (Vf(?!)i(rj), Vf02(ri)) " " " > Vf'?^iV£,(riV£,))^- is the dimension of D. 

The next quantity to evaluate is Some straightforward calculations yields: 

91n* _ 1 dD^ 1 dD^ dJ 

dCm dCm dCm 

The derivative of a Slater determinant respect to equals a sum of No determinants. 
The ith term of this sum is the determinant D with its iih. column Cj replaced with column 

y dCm ' dCm ' ' dCm > ' 

Eventually in order to evaluate |^ we proceed as follows (starting with (16)): 

fL^.W'-^-^^.^^^-'^ (21) 

The terms containing derivatives respect to can be evaluated as follows: 



dCm [D^^Y dCm dCm 

' + (22) 



{D^)^dCm dCm dc 

The second term gives the following expression: 

aVfln^ _ 2 dP^-^ 2 g^.pt . 

dc^ {D^rdcJ''' I (i)t)2 dc^ 

1 ^2^t+ 1 
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2 ^^V^^.|. 2 



•V^i^^ + 7^^^ + ^ (23) 



{Dl)^dCm ' dCm dc, 

1 (22) and (23) to evaluate and ^ 
W^D and and then implement the derivatives respect to 



Note that in equations (22) and (23) to evaluate ^J^^ and , we should first evaluate 
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